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Abstract 



We calculate the electron shear viscosity (determined by Coulomb electron collisions) for 
^ , a dense matter in a wide range of parameters typical for white dwarf cores and neutron star 

Q ' crusts. In the density range from ~ 10 g cm' -3 to 10^ - 10^° g cm ^ we consider the matter 

. composed of widely abundant astrophysical elements, from H to Fe. For higher densities, 

g cm ^, we employ the ground-state nuclear composition, taking into account 
finite sizes of atomic nuclei and the distribution of proton charge over the nucleus. Numerical 
values of the viscosity are approximated by an analytic expression convenient for applications. 
, Using the approximation of plane-parallel layer we study eigenfrequencies, eigenmodes and 

I viscous damping times of oscillations of high multipolarity, / ~ 500 — 1000, localized in the 

■ outer crust of a neutron star. For instance, at / ~ 500 oscillations have frequencies / ^ 40 

kHz and are localized not deeper than ~300 m from the surface. When the crust temperature 
decreases from 10^ K to 10^ K, the dissipation time of these oscillations (with a few radial 
^ ■ nodes) decreases from ~ 1 year to 10-15 days. 

1 Introduction 

o 

[ The shear viscosity of dense stellar matter (at densities p < 10^^ g cm~^) is important for a 
• number of astrophysical problems, including the viscous damping of oscillations in white dwarfs 
^ ' and in envelopes (crusts) of neutron stars. The total shear viscosity can be presented as a sum of 
^ . various matter components. In the outer crust of a neutron star or in the core of a white dwarf, it 
is determined by electrons and ions, rj = 7]^ + rj^. In the inner crust of a neutron star it is necessary 
to add a contribution of free neutrons, Under formulated conditions, the electrons are strongly 
degenerate and form an ideal Fermi gas, while the ions are fully or partially ionized and form a 
strongly nonideal Coulomb fluid or Coulomb crystal. In this case the electrons become the most 
important carriers of heat, charge (see, e.g., Ref. l^), and momentum. The main process of electron 
scattering which determines kinetic coefficients (thermal conductivity, electrical conductivity, and 
viscosity) is the scattering of electrons by ions (atomic nuclei). 

The shear viscosity of dense stellar matter determined by electron-ion scattering has been 
considered in a number of papers. For example, the electron viscosity of a strongly nonideal 
Coulomb liquid was calculated in Refs. [2 El El from variational principle. The results of these 
calculations are in good agreement. However, they were carried out neglecting the quasi-ordering 
of ions, which is important near the melting point. The inclusion of this quasi-ordering in a liquid 
together with multiple-phonon process of electron scattering in a crystal led to the disappearance 
of appreciable (by a factor of two to four) jumps in the electrical and thermal conductivities p]. 

Previous calculations of the viscosity were carried out in the Born approximation. However, 
the non-Born corrections are important for calculating the electrical and thermal conductivities 
of the matter containing chemical elements with high charge numbers Z (see, e.g., Ref. [T]). We 
include these corrections and show that they are equally important for calculations of the viscosity. 



While studying oscillations of neutron star envelopes, it is necessary to know the viscosity of 
the matter with the density p < 10^^ g cm~^. For p > 10^^ g cm~^, the sizes of atomic nuclei 
become comparable to the distances between them, and it is necessary to take into account the 
distribution of proton charge within the nuclei. This effect was included in the calculation of 
the electrical and thermal conductivity E] by introducing form factors of atomic nuclei. Such 
calculations have not been carried out yet for the viscosity. 

In this paper, we have calculated the shear viscosity taking into account non-Born corrections 
and the form factor of the nuclei, the quasi-ordering in a Coulomb liquid, and multi-phonon 
processes in a Coulomb crystal. The results are approximated by analytic formulas convenient for 
applications. 

Various types of oscillation modes can be excited in neutron stars. Generally speaking, these 
oscillations carry important information on the internal structure of neutron stars. Specific types 
of oscillations (such as r modes) can be accompanied by the radiation of gravitational waves. 
Since neutron stars are relativistic objects, their oscillations should be analyzed in the framework 
of General Relativity. The relativistic theory of oscillations was developed in a series of papers by 
Thorne and coauthors [3 El 13 UHl IHl 1121 • particular, the rapid (~ 1 s) damping of p-modes 
with multipolarity / = 2 owing to gravitational radiation was demonstrated in Ref. [2J. An exact 
inclusion of general-relativistic effects is complicated, but, in many cases, it is possible to use the 
relativistic Cowling approximation [IB]. A comprehensive analysis of various oscillation modes 
and mechanisms of their dissipation was carried out in Ref. . We also note the recent review of 
Stergioulas |15J, which contains an extensive bibliography. As a rule, oscillations with low values 
of / have been considered in the literature. 

Although neutron stars are in the final stage of stellar evolution, they can be seismically 
active for many reasons. Possible mechanisms for the generation of oscillations have been widely 
discussed in the literature (see, e.g., [^E] and references therein). Much attention has recently 
been paid to r modes — vortex oscillations that can be generated in rapidly rotating neutron stars 
and are accompanied by powerful gravitational radiation. In addition, oscillations can be excited 
in neutron stars, for example, during X-ray bursts (nuclear explosions on the surfaces of accreting 
neutron stars), bursting activity of magnetars (anomalous X-ray pulsars and soft gamma-ray 
repeaters; see, e.g., Ref. [in]), and glitches (sudden changes of spin periods) of ordinary pulsars. 

In this paper, we study the damping of oscillations in order to illustrate the viscosity calcula- 
tions. We choose the simplest example — p mode oscillations localized in the outer crust due to 
a high multipolarity, I > 500. 

2 Shear viscosity of dense stellar matter 
2.1 Parameters of dense stellar matter 

Let us consider dense stellar matter which consists of ions (atomic nuclei) and degenerate electrons. 
The state of strongly degenerate electrons can be characterized by their Fermi momentum p-p or 
wave number fcp, 

p-p = hk-p = h (Svr^ne)^''^ = meCXr, (1) 

where h is the Planck constant, rrie and are the mass and number density of electrons, Xr ~ 
1.009 (pqZ/A)^^^ is the relativistic parameter of the electrons, Z and A are the charge number and 
atomic number of the ions, and pe is the density in units of 10^ g cm"^. The electron degeneracy 
temperature is 

Tp = {ep - meC^) /^B ~ 5.93 ■ 10^ (^y^lTx^ - l) K, (2) 



where /cb is the Boltzmann constant and 



Ep = mlc =meCA/l + x^ (3) 

is the electron Fermi energy. In our study, we consider the matter at T ^ Tp and T < 5 ■ 10^ K 
(the latter is required to avoid dissociation of atomic nuclei). 
Let us also introduce the Fermi velocity of the electrons, 

vf = cPr = pp/ml. (4) 

The electrostatic screening of a test charge by degenerate electrons is described by the Thomas- 
Fermi wave- number k^F (the inverse electron screening length), 

where /x ~ is the electron chemical potential and a = e^/hc ~ 1/137.036 is the fine-structure 
constant. 

The state of ions is characterized by the Coulomb coupling parameter 

where a = [3/ (47mi)]^^'^ is the ion sphere radius; rij = rie/Z is the number density of ions; and Tg 
is the temperature in units of 10^ K. When F -C 1, the ions form a nearly ideal Boltzmann gas. If 
F > 1, they form a strongly nonideal Coulomb fluid. Finally, when F > Fm (which is realized at 
T < Tm), the ions crystallize. The crystallization of a classical system of ions occurs at Fm ~ 175 
(see, e.g., [E]). 

Quantum effects in ion motion become important at = T/Tp <^ 1, where 

Tp = hujp/kB ^ 7.832 ■ 10%Z/A)pI^'^ K (7) 

1 /2 

is the ion plasma temperature, Up = {ATxZ'^e'^n\/m\) is the ion plasma frequency, m\ = Am^ is 
the ion mass, and mu = 1.6605 x 10~^^ g is the atomic mass unit. 
For isotropic matter, the viscous stress tensor has the simple form 

, _ fdU„ , dUi, 2 



where U is the hydrodynamical velocity of the matter, rj is the shear viscosity, and ( is the bulk 
viscosity (this last quantity is especially important for uniform compressions and rarefactions of 
the matter). 

Generally, crystalline matter is anisotropic, and Eq. (jHI) for the viscous stress tensor may 
be formally invalid. However, in a dense matter, ions crystallize into highly symmetric body- 
centered cubic (bcc) or face-centered cubic (fee) lattice. In this case, the viscous stress tensor 
for a monocrystal is determined by three independent coefficients (see, e.g., Ref. |IH|); and can 
be written in the form (jHJ with an additional term of the form KdapdUa/dxa (the sum over a 
is not carried out). The quantity U should be treated as the velocity field for shifts of the ions 
in their lattice sites. When studying transport processes on scales exceeding the characteristic 
monocrystal size, the matter can be considered as isotropic. As in all the literature concerned 
with the kinetics of crystalline matter of white dwarfs and neutron stars without magnetic fields, 
we will restrict our analysis to this case (assuming k = 0). 



The shear viscosity of neutron star crusts and white dwarf cores is primarily determined by 
strongly degenerate electrons. It is convenient to present this viscosity in the form 



* = (9) 

5z/e 

where i^e = l/^e is the effective electron collision frequency, which is the inverse of the effective 
electron relaxation time t^- If the electron scattering is determined by several independent pro- 
cesses, these can be studied separately, and the total collision frequency will be the sum of partial 
ones. For the dense matter of white dwarf cores and neutron stars envelopes, one has 

= ^ei + J^imp + ^ee, (10) 

where z/gj, t'imp and z/ge correspond to electron scattering by ions, impurity ions, and electrons, 
respectively. The dominant process is electron-ion scattering; it will be studied below. This 
scattering determines also the thermal and electrical conductivities of dense matter (see, e.g., 
Ref. [T]). With small variations, the formalism proposed by Potekhin et al. is applicable for 
calculations of the viscosity. 

In crystalline matter, the electron-ion interaction can adequately be described in terms of the 
emission and absorption of phonons ^H] • This description can be realized using an ion dynamical 
structure factor 0. 

The frequency of electron-ion collisions {ei scattering) can be written as 

I/gi = 127r 2 = ^ei, (11) 

where Agi is the effective Coulomb logarithm, which can be calculated using the variational method 
(see, e.g., Ref. [Hj). Using the simplest trial function in the Born approximation for a strongly 
nonideal ion plasma (F > 1), one obtains 
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5,(g)dg, (12) 



where go is the minimum momentum q transferred in an ei scattering event; go = for the 
liquid phase and go = ge in the crystalline phase, where ge = (Ovrni)^^'^ is the radius of a sphere 
of the same volume as the Brillouin zone. The value go = ge selects Umklapp processes (i.e., 
those involving electron momentum transfer g ^ ge) in an ei scattering event. At not very low 
temperatures 

T>T„~TpZi/3a^/(3/3,), (13) 

the contribution of such processes to the Coulomb logarithm is much higher than the contribution 
of normal processes which occur at g < ge (see, e.g., Ref. [201 )• However, at low temperatures 
(T < Tu), Umklapp processes are "frozen out" and the viscosity is determined by normal processes. 
We will neglect this effect below, restricting our consideration to temperatures T > T^. 

The function u{q) in Eq. ()12|) describes the Coulomb interaction between an electron and 
an atomic nucleus; it is discussed in Section 2.3. The factor in the square brackets describes the 
kinematic effect of backward scattering of relativistic electrons (see, e.g., |2I]); ^^(g) is an effective 
ion static structure factor that takes into account ion correlations. This factor coincides with the 
effective structure factor that determines the electrical resistivity of the dense matter (it was 
calculated and approximated in Ref. |22])- Note that the structure factor of a strongly nonideal 
Coulomb liquid is known only in the classical limit (B ^ 1). We also define a simplified structure 
factor, based on the two approximations: (1) Neglecting quasi-ordering in ion positions in the 



Coulomb fluid (see, e.g., Ref. |T]). (2) Single-phonon approximation for the inelastic structure 
factor of the Coulomb crystal (see, e.g., Ref. [2111 )• We will call the viscosity calculated using the 
simplified structure factor the simplified viscosity. Note that previous calculations of the shear 
viscosity by Flowers and Itoh |2j sBj and by Nandkumar and Pethick [31 were carried out for a 
Coulomb liquid using the simplified structure factor. 

To account for non-Born corrections, we also multiply the integrand by the ratio of the exact 
and Born cross sections for Coulomb scattering. This method was proposed in Ref. and used 
for calculating the transport coefficients by Potekhin et al. IT]. 

The effective frequency of electron scattering by impurities (assuming that impurity atoms are 
randomly distributed over the crystal) is similar to the frequency of electron-ion scattering (see 
Eq. (HH)): 

^imp 9 / , -^imp) ^imp ^impi (-^^) 

^ imp 

where Z\^^ is the charge number of an impurity ion and the Coulomb logarithm Ajmp is calculated 
using Eq. (fT^ . but assuming that impurities are weakly correlated (which corresponds to the 
structure factor Simp = 1, while the screening of impurities is included into the factor u{q)). In 
the simplest model of Debye screening (with a screening length q^-^p) 

Aimp = ^ (1 + + 2^ + 2^(3^) In (l±i) -^-(3^^-^-f3^- 1, (15) 

where ^ = qsimp/ {2kp); qsi^ip = ^tf + ^fmp? ^^^1 fcimp is the wave number for the Debye screening 
of a test charge by impurities (the inverse correlation length of impurities). This wave number 
weakly affects the Coulomb logarithm (fcxF ^ ^imp), and can be estimated as k^^p = (47minip/3)^/^, 
where n^^p is the number density of impurities. Electron-impurity scattering is important at low 
temperatures, when electron-phonon scattering is suppressed by quantum effects. 

The expression for the frequency of electron-electron collisions z/ge was obtained by Flowers 
and Itoh Their result can be written in the form 



^■K^a^klT"^ f kF \ 6 2 \ , 1 nii W ^0 



9*2, V, / I- . 2 + ^)- 4.473.10^M-^) (::H)'''t| S-, (16) 

where the latter expression is presented for an ultra- relativistic electron gas (xr ^ 1), no ~ 
0.16 fermi"^ is the number density of nucleons in atomic nuclei, and Tg is the temperature in units 
of 10« K. 

2.2 The form factor of atomic nuclei 

The function u{q), which describes the Coulomb electron-ion interaction in Eq. (fT^ has the form 

q^\e{q)\ 

where e{q) is the static longitudinal dielectric function of the degenerate electron gas [23], which 
takes into account the electron screening of the ion charge. Here, 

F{q) = ^ j enp{r) exp(^ r-q)dV=^ j ^p(^) ^^^(g^) ^ dr (18) 



is the nuclear form factor that characterizes the distribution of proton charge within an atomic 
nucleus. The integration in Eq. (fTH]) is carried out over the atomic nucleus, np{r) is the local 



number density of protons, and is the radius of the proton core. In white dwarfs and the outer 
envelopes of neutrons stars (p < 10^^ g cm~^), the atomic nuclei can be considered as pointlike, 
with F{q) = 1. At densities p < 10^^ g cm~^, the proton charge can (with good accuracy) be 
taken uniformly distributed throughout the nucleus. In this case, the form factor is 

3 

Fiq) = -^^^ [sin(grp) - grpcos(grp)] 

For p > 10^^ g cm~^, the proton density profile strongly differs from a step function, and the form 
factor (fT^ becomes unacceptable. Then we have employed the nuclear form factor provided by the 
model of the ground-state matter (with the parameters smoothed over density jumps associated 
with the changes of nuclides [2S1)- 



(19) 



2.3 Analytic approximation for the viscosity 

We have obtained the analytic approximation for the Coulomb logarithm of ei scattering using the 
method of the effective electron-ion scattering potential proposed in Ref. [Tj for the electrical and 
thermal conductivities. The conductivities of the matter at p < 10^° g cm~^ were studied in P, 
where the form factor of atomic nuclei was taken to be unity. Later, Gnedin et al. ^27] extended 
this method to higher densities using an approximate nuclear form factor. As mentioned above, 
the effective structure factors for the viscosity and electrical conductivity coincide. This simplifies 
the generalization of the effective-potential method for describing the shear viscosity. Following 
Ref. [21!, we replace u^{q)S^{q) in Eq. ^ by 
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(20) 



The factor (g^ + 93) ^ corresponds to the Debye screening of the Coulomb interaction with the 
effective screening length gg^; the term in the square brackets describes ion correlations. The 
functions and D characterize ion quantum effects. The factor exp (— Sig^) added in Ref. j2ZI 
takes into account the effect of the nuclear form factor. Numerical values of the shear viscosity 
obtained from the exact theory are reproduced taking the same parameters as for the electrical 
and thermal conductivities in Ref. [77j . 
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where Se = k'^p/{2kF)'^ = a/{TT(3y) is the electron screening parameter, ro = a/\/^ is the ion 
Debye screening length, x^^c is the ratio of the mean-square radius of the proton distribution 
within an atomic nucleus to the ion sphere radius; and m_i ~ 2.8 and m_2 ~ 13 are the frequency 
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Figure 1: Temperature dependence of the shear viscosity for a carbon plasma with the density 
p = 10"^ g cm~'^ . The sohd curve is the analytic approximation. Bold points present numeri- 
cal calculations. The dashed curve is the "simplified" viscosity, which demonstrates a jump at 
the melting point. The dot-dashed curves correspond to the matter with ^^O impurities with 
concentration of 2 and 4%. 



moments of the phonon spectrum in the bcc Coulomb crystal. Note that the function coincides 
with the function in Ref. p71 . 

After integrating in Eq. (|12|) with the effective potential (|2m) . we obtain 

A=[Ao{s,w + wi)-Ao{s,wi)]Gr,D, (27) 

where 

Ao(s, w) = Ai(s, w)-{l + p^) A2{s, w) + P^Asis, w), (28) 



2Ai{s,w) = + [l-e-"") -{1 + 310)6'"" [Ei{sw)-Ei{sw + w)] 

2A,is,w) = ^""-^ + ^ -^(l-e-)-2. In^ 
w s + 1 s 

+s{2 + s w) e'"" [El (sw) - Ei {sw + w)] , 

2 , 1 + s 1 2s3 _ - 3s + 1 ,„ e""' 

2A3S,u; = 3sMn + - re-"' + 



2 1 + s (l + s) w 

(2sw- l)(l-e "-) _ ^2 ^3 ^ ^Ei{sw) - Ei{sw + w)) . 



In this case, Ei{x) = J°° y ^e ^ dy is the exponential integral (see, e.g., Ref. |2H1)- The maximum 
error in our approximation of the viscosity does not exceed 20%. 



2.4 The main properties of the shear viscosity 

Let us discuss the results of our calculations of the shear viscosity neglecting the freezing out 
of Umklapp processes (Section 2.2). Figure Q presents the temperature dependence of the shear 
viscosity for a carbon plasma at p = 10^ g cm~^. The upper horizontal axis shows the Coulomb 
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Figure 2: Temperature dependence of the shear viscosity for an iron plasma with the density 
p = 10^ g cm~^ . The sohd curve is the analytic approximation. Bold points present numerical 
calculations. Hollow circles correspond to the Born approximation. The dashed curve is the 
"simplified" viscosity calculated in the Born approximation. 

coupling parameter F of ions. Because the charge number is rather low, Z = 6, non-Born cor- 
rections are small and invisible in Fig. ^ All the data presented in the figure, except for the 
dot-dashed curves, correspond to scattering of electrons by ions of single type. 

Bold points in the figure show numerical results. The solid curve is the analytic approximation. 
The dashed curve is the viscosity calculated using the simplified structure factor (Section 2.2). 
Large jumps of this "simplified" viscosity at the melting point are clearly visible. These jumps 
(by a factor of two to four) are present for all chemical elements and all plasma parameters. Our 
modification of the structure factor (Section 2.2) increases the accuracy of calculations in the 
liquid and solid phases, and nearly removes the viscosity jumps for all elements. This makes it 
possible to introduce a unified approximation for both phases (Section 2.4). 

Nevertheless, appreciable viscosity jumps are still present at the melting point in our calcula- 
tions for high densities, where zero point oscillations of ions become important. We assume, as 
did Potekhin et al. [1^ for the electrical and thermal conductivities, that these jumps result from 
using the classical structure factor in the ion liquid under the conditions when quantum effects 
are important (because quantum effects are correctly included only in the solid phase). Since 
numerical data used to construct the analytic approximation include both the liquid and solid 
phases, the unified analytic approximation shifts the viscosity in the liquid phase to the viscosity 
in the solid phase. We expect that, for an ion fluid at high densities, this approximation is more 
exact than our original numerical data. It will be possible to verify this in the future, when ion 
structure factors in a fluid are calculated taking into account quantum effects. 

The dot-dashed curves in Fig. [T] demonstrate the effect of scattering by charged impurities. 
We have considered oxygen impurities with concentrations of 2 and 4%. The presence of these 
impurities is weakly manifested at high temperatures, but dominates at low temperatures, T <^ Tp, 
when electron-phonon scattering in the Coulomb crystal is strongly suppressed by quantum effects. 

Figure |21 shows the temperature dependence of the shear viscosity for an iron plasma at p = 
10^ g cm~'^ . The upper horizontal scale plots the Coulomb coupling parameter F. Bold points 
present our numerical results, while the solid curve shows the analytic approximation. The dashed 
curve is calculated with the simplified structure factor neglecting non-Born corrections. As in 




Figure 3: Shear viscosity at the temperature T = 10'' K as a function of the density for various 
chemical compositions (H, He, C, and Fe). Bold points are the numerical results and the curves 
show the analytic fit. 

Fig. m the simplified viscosity displays jumps at the melting point, while the new results pass 
smoothly through this point. The charge number of iron {Z = 26) is high enough for the non- 
Born corrections to be appreciable. To demonstrate this effect, hollow circles in Fig. |21 show the 
calculated viscosity in the Born approximation. We can see that the non-Born corrections reduce 
the viscosity by approximately 20%. 

Figure El shows the density dependence of the shear viscosity for hydrogen, helium, carbon, 
and iron plasmas at T = 10'' K. We consider the densities typical for the cores of white dwarfs 
and the outer envelopes of neutron stars. Bold points show numerical results, and the curves are 
analytic approximations. The strong dependence of the viscosity on the chemical composition is 
due to the dependence of an electron-ion collision frequency on the charge number Z. In contrast 
to the thermal conductivity (see, e.g., Ref. [Tj), the effect of electron-electron colhsions on the 
viscosity is insignificant at the considered densities, even for hydrogen. 

Figure m demonstrates the density dependence of the shear viscosity in the density range from 
10^ to 10^^ g cm~^for the three temperatures T = 10'', 10^, 10^ K. The ground-state nuclear 
composition is employed with smoothed parameters. Symbols show the numerical results, and 
curves are analytic approximations. In contrast to the thermal conductivity jT] , the shear viscosity 
decreases strongly with growing temperature. Note that the ratio p/rj grows with increasing 
density in the outer crust of a neutron star. These results are important for the damping of 
oscillations in the neutron star crust (see Section 3.3). 

For illustrative purposes, the plot is continued beyond the crust into the stellar core (to densities 
p > 1.5 ■ 10^^ g cm~^ ). In the core, we have used the equation of state of the matter presented in 
Ref. ^j. It is assumed that the core is composed of neutrons, protons, and electrons and is not 
superfluid. The electron viscosity in the core of such a star is mainly determined by the scattering 
of electrons by degenerate protons. The corresponding collision frequency is obtained in the same 
way as the electron-electron collision frequency (see Eq. (fTH|) ). 

..P = .V ('-l) 1.434. 10- f^E) TSC-^X '-^ c-\ (29) 



where rrip is the proton mass, and m* is its effective mass, which differs from rrip due to manybody 
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Figure 4: Shear viscosity of the ground-state matter versus p for the three temperatures T = 10'', 
10^, 10^ K. Sohd curves show the analytic approximation; symbols are our numerical calculations. 
The vertical dotted lines indicate the neutron drip density and the crust-core interface in a neutron 
star (set to be p = 1.5 x 10^"^ g cm~^ ). The electron shear viscosity in the stellar core, determined 
by scattering of electrons by degenerate protons, is presented for comparison. 

effects (we have set m* = 0.7 mp). The Debye-screening parameter in the stellar core is 



where the sum is taken over all types of charged particles (electrons and protons); Cj, Uj and /ij 
are the charge, number density, and chemical potential of particle species j. Owing to a strong 
suppression of the scattering rate by proton degeneracy, the electron viscosity in the core grows 
by a factor of ~ 1000 as compared to the viscosity in the crust. 

3 P mode oscillations of a neutron star crust 

This section is devoted to p mode oscillations (i.e., oscillations in which perturbations of the 
pressure dominate over the buoyant force) with high orbital numbers (multipolarity), I > 500, 
localized in the outer crust of a nonrotating neutron star. 

3.1 General formalism 

3.1.1 Flat metric for the envelope of a nonrotating neutron star 

The spacetime metric for a nonrotating neutron star 30 can be written as 



where df2^ = d6'^ -|- sin^ 6'dv9^, t is the time coordinate, r is the radial coordinate, d and are 
the polar and azimuthal angles, and the functions A(r) and $(r) determine spacetime curvature. 
It is sufficient to consider a thin envelope, neglecting the variation of the functions A(r) and $(r) 




(30) 



ds2 = cV*d?-e2"dr2-r2dfi2, 



(31) 



in the crust and using their values at the stellar surface, 



where M is the gravitational stellar mass. Neglecting variations of r in the envelope compared to 
the stellar radius R (in the approximation of a thin envelope layer), we can rewrite Eq. (|31|) in 
the form ^ 

= c'(^l- dp -(^1-^^ dr' - dn\ (33) 

where Rq = 2GM/c^ ^ 2.953 (M/Mq) km is the gravitational radius. Introducing the local time 
t and the local depth z, 

t = i^/l-RG/R, z={R-r)/^l-RG/R, (34) 
we come to a flat coordinate system in the neutron star crust, 

ds^ = c^dP- dz^-R^dn\ (35) 

3.1.2 Equilibrium structure of the crust 

The structure of a neutron star is determined by the equation of hydrostatic equilibrium, including 
the effects of General Relativity (the Tolman-Oppenheimer-Volkov equation; see, e.g., Ref. [5U]). 
This equation is greatly simplified in the envelope; in the coordinate system (jMjl it can be written 
as 

where Pq and po are the equihbrium pressure and density, = dPo/dpo is the square of the local 
sound speed, and 



GM ^ M /lOkmA^ / cm 



9 = r^'^ - 1-327- 10^^-^ I ) /JI-Rg/R ^ (37) 

R^^I-Rg/R Mq \ R ^ ^ ' 



is the gravitational acceleration. 

In computations we have used the equation of state of a fully degenerate electron gas with 
electrostatic corrections and have employed the model of the ground-state matter with smoothed 
parameters. Alternatively, we have used a polytropic envelope, composed of ^^Fe nuclei, where 
the pressure is determined by degenerate electrons, assumed to be relativistic at all densities. 

3.1.3 The oscillation equation 

In the outer envelopes of neutrons stars, the main contribution to the pressure is produced by 
degenerate electrons. Therefore, while considering p-modes, we can use a single equation of state 
to describe the equilibrium configuration of the star and perturbations. 
Let us use the Euler equation in a fiat-space metric (jH^j) . 

dU VP 

+ {U.V)U = -^+g, (38) 
ot p 

where P is the pressure. The continuity equation must also be satisfied, 

^ + V (pU) = 0. (39) 



Assuming that the velocity U is small and introducing Euler perturbations of the pressure SP = 
P — Pq and the density 6p = p — po, we obtain the linearized Euler equation 



^ = ^fvPo--V«P, (40) 

dt pi Po 

and the continuity equation 

^ + V (poC/) = 0, (41) 
while the equation of state for the perturbations can be rewritten in the form 

5P = cl5p. (42) 

We will consider irrotational motion and write the velocity in the form U = V0, where is the 
velocity potential, which is a scalar function of coordinates and time. Formally, the function 
is determined up to an arbitrary function of time, which we choose in such a way that the Euler 
equation can be rewritten (using Eqs. (jHBj) and ((121)) as 

^ = = -c^^ (43) 
dt Po ^ ^ 

Differentiating ()43p with respect to t and taking into account Eqs. (PT|) and ()36|) . we get 



ciA<j) + g-V<j), (44) 
where we have introduced the Laplace operator 

dz^ ^ [de^ ^ sin^ edcpy ' ^ ' 

An equation that coincides with (j44j) was obtained by Lamb [21] for oscillations in the Earth 
atmosphere. The variables in (|44j) can be separated if we write 

(j) = e'^'YiUO,v)F{z), (46) 

where iv is the oscillation frequency, Yim{0, ip) is a spherical harmonic function (see, e.g., Ref. [S2]), 
and F{z) is an unknown function of the depth determined by 

'^".^^.(^-^V-O. (47) 



dz"^ Cg dz \ Cg , 

The first boundary condition for this equation, 

function F{z) is bounded as z — > 0, (48) 

follows from the requirement of a finite oscillation amplitude at the stellar surface. The second 
boundary condition is imposed artificially. In the current study, we have to solve equations that are 
applicable only in the thin stellar crust. Therefore, oscillations should be damped at sufficiently 
high depths. For simplicity, we formally move this boundary condition to infinity along z and 
will check a true localization of oscillations (see Section 3.3). In this case, the second boundary 
condition can be written 

F{z) ^0 as z ^ oo. (49) 



Together with the boundary conditions (jiHjl and (jl^ . Eq. (jTTj) specifies eigenfrequencies and 
eigenmodes of oscillations. Moreover, the following asymptotes are valid at large and small depths: 



Fi^)^{ ..^^f "mTT^.fA (50) 




+ 1) z/R^ z 



oo. 



Perturbations of the pressure and the density are expressed in terms of the function 0(r) using 
Eq. iSD, 

6P = -too po(j), 5p = — = -i^^ (j). (51) 



Due to the boundary condition (jiH|) . variations of the pressure and the density, 5P and 5p, vanish 
at the stellar surface (because po{R) = 0). We can see from these last expressions that the number 
of radial nodes {k) of the velocity potential coincides with the number of nodes of the pressure 
and density variations. Below, we will call k the number of radial nodes of the mode. 
The displacement vector for an oscillating matter element can be written in the form 



The z component of this vector is 



Udt = ■ (j). (52) 

UJ 



1 dF 

i. = ~-Yi^{e,^)—, (53) 
UJ az 



and the magnitude of horizontal displacement can be estimated as 

l^hl^-^ \F{z)\. (54) 

UJ It 

The quantities lF{z)/R and dF/dz appear in Eq. ()47|) on equal footing. Therefore, horizontal 
and radial displacements should have the same order of magnitude for the oscillations considered. 

Oscillations of a polytropic envelope in a plane-parallel approximation were studied earlier by 
Gough jSS] assuming that the equations of state of the unperturbed matter and the perturbations 
are polytropic with different indices. In the limiting case of the same polytropic index n, his result 
can be presented as follows. The mode containing k radial nodes has the eigenfrequency 

0^1 = I (I + l) « 10»9u (^) VWTY, (I + l) (55) 

while the velocity potential is given by 

Fk{z) = exp {-VWTT)^) Lf^ (2 V7(7Tl) I) , (56) 

where L^" ^\x) is a generalized Laguerre polynomial (see, e.g., Ref. |2H]) and (?i4 is the gravita- 
tional acceleration at the stellar surface in units of 10^'* cm s~^. Note that the eigenfrequencies 
agree with the simple estimate u"^ ~ g/a, where a ~ R/l is the characteristic depth for the 
oscillation localization. 

Let us remark that the mode with k = 0, which does not have any radial nodes, corresponds to 
the vanishing Lagrangian variations of the pressure and the density (V- U=0, see (f7T|) ): its param- 
eters do not depend on the polytropic index. Adding the condition V ■ C/ = A0 = to Eq. 

it is easy to show that this mode is described by the function F{z) = exp (^-y/l{l + l)z/Rj and 

has the frequency 

^0' = I VWTT) ^ 10« g^, (i^) VWTT) s-^; (57) 



it exists for any equation of state. Therefore, the frequency c^o will further be used to normalize 
oscillation frequencies. 

The frequencies uj we discuss refer to the coordinate system of the stellar envelope (see 
Eq. ()34j)). They are easily transformed to the frequencies uo detected by a distant observer, 

u = uj Vl - Rg/R- (58) 



3.2 Viscous damping of oscillations 

In this section, we consider the damping of oscillations with velocity potentials of the form 
e*"^* ^m(6', v^)-F(-2) in a spherically symmetric star under the action of shear viscosity. We take 
flat space-time metric which is applicable to oscillations studied in Section 3.1. This is because 
the flat metric which coincides with the metric for a thin spherical layer in a fiat space-time, 
can be introduced in the region, where oscillations are localized. As a result, it is sufficient to 
consider the oscillation damping time in a flat metric and transform this time (in accordance with 
Eq. into the frame of a distant observer. We define the oscillation damping time r as 



where 



r = E/\dE/dtl 



E= I £dV = I p^-^dV 



is the total energy of oscillations and 



^ = Upo\U\' + ^\6p\' 
4 V Po 



(59) 
(60) 

(61) 



is the energy density of oscillations at a given point averaged over an oscillation period (see, e.g., 
Ref. ^j). The additional factor of 1/2 in the expression for e is required due to the averaging 
over the oscillation period. The integration is carried out over the entire volume of the star (in 
practice, over the region where the oscillations are localized). We neglect perturbations of the 
gravitational potential. The last equality in Eq. (jUUI) is determined by the equality of the mean 
kinetic and potential energies for low-amplitude harmonic oscillations. Note that some authors 
introduce the damping time for the oscillation amplitude rather than the damping time (j^^ for 
the oscillation energy. 

When calculating the energy from Eq. (j6U|) . the angular integration can be carried out analyt- 
ically. 



E 



R 



P 



(F') 



r dr. 



The period-averaged rate of the viscous energy dissipation is (see, e.g., Ref. f^) 



dE 
~dt 



du* dui'. 

" ' dK 



+ 



dxi 



(62) 



(63) 



where the viscous stress tensor cr^^ is given by Eq. (jH}. As in the expression for the oscillation 
energy E, the additional factor of 1/2 is required owing to averaging over the oscillation period. 
It is easy to see that the energy dissipation rate separates into a sum of terms associated with 
shear and bulk viscosities. We will restrict ourselves to the dissipation determined by the shear 
viscosity. The integration over the angular variables in Eq. (jU^ can be carried out analytically 
(see the Appendix). 
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Figure 5: Eigenfrequencies of oscillations localized in the crust of a "canonical" neutron star. The 
frequencies are normalized to the frequency ujq given by Eq. (jHTjl . The numbers next to the curves 
indicate the number of radial nodes. 

3.3 Discussion of numerical results 

As an example, let us choose a "canonical" model for a neutron star with the mass M = IAMq 
and the radius R = 10 km. For this model, 



The thickness of the outer crust in such a star (p < 4 ■ 10^^ g-cm^'^, before the neutron drip point) 
is ~ 364 m. 

Oscillation eigenfrequencies have been found via a series of iterative trials, checking for the 
coincidence of a mode number with a number of radial nodes. 

3.3.1 Oscillation eigenfrequencies 

The dependence of oscillation eigenfrequencies on /, calculated from Eq. (jlTj) with the boundary 
conditions ()48j) and (jlHI), is presented in Figs. El and El As mentioned earlier, the frequency of 
the fundamental mode, which does not have any radial nodes, is determined by Eq. ()57|1 for all /. 
With decreasing /, oscillations penetrate deeper into the outer crust, where the equation of state is 
softened because the electron gas becomes relativistic and electrons undergo beta captures. This 
slightly decreases the dimensionless eigenfrequencies. As in the model with the polytropic equation 
of state (see Eq. (jHSl)), the separation between squares of neighboring eigenfrequencies for any / 
is nearly constant. A weak decrease of this separation with the growing number of radial nodes 
is due to the penetration of oscillations into deeper layers of the star, where the equation of state 
is softer. When I ~ 500, the main oscillation energy is localized in the region 50 m< z <400 m, 
where the equation of state of the degenerate, relativistic electron gas well describes the profile of 
the sound speed. Therefore, the estimate (j55|l for eigenfrequencies of the polytropic crust model 
with the polytropic index n = 3 is justified (deviations constitute several per cent). 




(64) 



and for a distant observer 




(65) 
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Figure 6: Frequencies of oscillations localized in the crust of a "canonical" neutron star as detected 
by a distant observer. The subscript of / denotes the number of radial nodes. 
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Figure 7: Root-mean-square (over angles) amplitude of radial displacements of matter for modes 
with I — 500. The subscript of ^ indicates the number of radial nodes. The root-mean-square 
amplitude of radial displacements at the stellar surface has been set equal to 1 m. 



3.3.2 Oscillation eigenmodes 

Figure m presents profiles of radial displacements of matter elements for the modes with / = 500. 
The root-mean-square amplitude of radial displacements of the stellar surface has been set 1 m. 
Since we are studying linear oscillations, this quantity is an arbitrary (sufficiently small) constant 
that normalizes the solution. From Fig. [7| one can easily determine the magnitude of radial 
displacements in the star for any other amplitude of the displacements at the surface. With 
growing z, the radial displacement amplitude overall decreases. At z ~ 300 m, the decrease 
becomes monotonic and gradually tends to the exponential asymptotic (jKn|) . This means that 
oscillations are localized in the outer crust. 




z[m] 

Figure 8: Angle- averaged energy density e of oscillations with / = 500. The subscript of e indicates 
the number of radial nodes. The root-mean-square amplitude of radial displacements at the stellar 
surface is 1 m. 

This effect is manifested even more clearly in the energy density of oscillations. Figure |H1 
presents the dependence of the angle-averaged energy density e versus depth z for oscillations with 
/ = 500. The mode amplitude is normalized in the same way as in Fig. [3 In our approximation, 
the energy density is proportional to the square of the normalized amplitude of radial surface 
displacements. The depicted modes are localized in the outer crust of the star. The energy 
density of oscillations varies relatively weakly within the "critical" depth z < 100 — 200 m, after 
which it falls off exponentially. The energy density decreases by more than two orders of magnitude 
near the boundary between the outer and inner crusts. 

As noted in Section 3.3.1, at / ~ 500 oscillation frequencies are well reproduced by a polytropic 
crust model. The situation is somewhat different for oscillation eigenmodes. Normalization at the 
stellar surface is not expedient for these modes, since the polytropic model poorly reproduces 
the structure of the star at low depths, z < 40 m. Consequently, such a normalization leads to 
large errors at the depths of interest to us, z < 100 — 200 m, where the main oscillation energy 
is concentrated. Therefore, we need a special normalization to compare exact and polytropic 
modes. Figure El (similar to Fig. ^ shows the angle-averaged energy density e of oscillations 
as a function of z. The symbols refer to exact numerical results, while the curves are for the 
polytropic model, normalized in such a way to be in agreement with the exact results in the 
region where oscillations are localized. We can conclude that the polytropic model for the outer 
crust satisfactorily reproduces the energy density of oscillations at depths of 60 m < z < 300 m 
for modes with / < 500. 




z[m] 

Figure 9: Angle- averaged energy density e of oscillations with / = 500. The subscript of e indicates 
the number of radial nodes. The curves refer to the polytropic model, and the symbols are for the 
exact equation of state. 
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Figure 10: Damping time r of oscillations for a distant observer as a function of multipolarity / 
for a neutron star with the internal temperature T = 10^ K. The subscript of r gives the number 
of radial nodes. 



3.3.3 The oscillation damping 

In our calculations we have considered the neutron star crust as isothermal. This approximation 
describes well the internal temperature profile. The temperature is nearly independent of depth 
due to a high thermal conductivity of degenerate electrons. Naturally, oscillation frequencies and 
damping times do not depend on the normalization amplitude (the amplitude of radial surface 
displacements). Figures ITIH ITT] and IT^ present the dependence of the oscillation damping time r 
(for a distant observer) on / for a canonical neutron star with the internal temperature T = 10'', 
10®, and 10^ K. The strong temperature dependence of the damping time is due to the appreciable 
decrease of the viscosity with increasing temperature (Fig.EJ. The damping time can be estimated 
from characteristic oscillation parameters as 

r~£/£~pf/2 / hf^V) -A^^, (66) 




where e is the local viscous dissipation rate, while U and A are, respectively, the characteristic 
velocity and its variation length scale in the region of oscillation localization. 
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Figure 11: Same as Fig.^lbut for the internal temperature T = 10*^ K. 
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Figure 12: Same as Fig. [101 but for the internal temperature T = 10 K. 



Let us consider Fig. ^2 in more detail. Oscillations with / ~ 500 are localized at 2; < 100 m 
(Fig. IHI), which corresponds to densities p < 10^° g cm~^ . Under these conditions, the ratio p/rj is 
~ 3 s cm~^ (Fig. and grows with increasing / (due to the decrease of the density in the region of 
oscillation localization). Let us make estimates for the modes with I ~ 500. The velocity variation 
length scale can be estimated as A ~ R/l. Note that this scale decreases for modes with a large 
number of radial nodes, accelerating the damping of these oscillations. For the fundamental mode 
(without any radial nodes), the damping time (in the frame of a distant observer) can be estimated 



as 



r 



1.2 10^ X (500//)^ s ^ 120 X (500//)^ day. 



(67) 



in good agreement with the numerical results for I ~ 500. The damping time drops weaker than 
oc because the ratio p/r] grows for higher /, owing to the decrease of the density in the region 
of oscillation localization. 

For the outer crust of a star with the temperature T = 10^ K, the ratio p/rj depends weakly 
on p. Therefore, the oscillation damping time obeys the law r oc 

The damping time as function of the crust temperature is presented in Fig. where we have 
chosen the modes with / = 500 as an example. The damping time grows by approximately two 
orders of magnitude as the temperature varies from 10^ to 10^ K. When the temperature increases 
by another order of magnitude, the damping time grows further by a factor of three. This is due 
to the nonlinear decrease of the viscosity with the growth of T (Fig. 0}. 

There are many other oscillation damping mechanisms in addition to the viscous damping 




7,0 7,5 8,0 8,5 9,0 

ig 7-(K) 

Figure 13: Damping time of oscillations as a function of the internal neutron star temperature for 
modes with / = 500. The subscript indicates the number of radial nodes. Numerical results are 
redshifted for a distant observer. 

considered here. For example, one often considers the damping due to the emission of gravitational 
and electromagnetic waves (for oscillations of the stellar matter with a frozen-in magnetic field). 
In our case, these mechanisms are inefficient due to high multipolarity, / ~ 500. If / is high, 
we expect gravitational or electromagnetic radiation to be generated by an ensemble of closely 
spaced coherent elementary radiating regions, which radiate in antiphase and cancel one another. 
Formally, the weakness of this radiation is manifested by the presence of large factors {21 + 1)!! 
in the denominators of the expressions for the radiation intensities (see, e.g., Refs. An 
analysis shows that the damping of oscillations we consider here is determined to a great extent 
by the shear viscosity. 

A detailed analysis of the evolution of pulse shapes of some radio pulsars provides an evi- 
dence that high-multipole oscillations are, indeed, excited in them (see, e.g., the recent paper 
|37j). However, reliable observational data on the existence of such oscillations have not yet been 
obtained. 

4 Conclusions 

We have calculated the shear viscosity of dense stellar matter for a broad range of temperatures 
and densities typical for the cores of white dwarfs and the envelopes of neutron stars. We have 
considered the matter composed of astrophysically important elements, from H to Fe, at the 
densities from 10^ — 10^ g cm~'^ to 10^ — 10^° g cm^^ . At higher densities, 10^° — 10^^ g cm~^ , we 
used the model for the cold ground-state matter taking into account finite sizes of atomic nuclei 
and the distribution of proton charge within the nuclei. Under the conditions described above, 
the shear viscosity is determined by Coulomb scattering of degenerate electrons by atomic nuclei. 
We have used the modified structure factor of ions proposed by Baiko et al. j22] and applied by 
Potekhin et al. IJ for calculating the thermal and electrical conductivities of dense matter. In the 
ion liquid, this modification approximately takes into account a quasi-ordering in ions positions, 
which reduces the electron-ion scattering rate. In the crystalline phase, the new structure factor 
takes into account multi-phonon processes, which are important near the melting temperature T^. 
The new results near the melting point differ appreciably from the results of Flowers and Itoh 
12 EI obtained for a Coulomb liquid. We have approximated our numerical results by an analytic 
expression convenient for applications. 

We have analyzed eigenfrequencies and eigenmodes of oscillations localized in the outer crust 
of a neutron star (in the plane-parallel approximation). A polytropic model for the crust can 



reasonably well reproduce the eigenfrequencies of oscillation modes with multipolarity I ~ 500. We 
have also calculated viscous damping times of the oscillations. We have obtained a sharp decrease 
in the damping time with decreasing the internal temperature of a neutron star. For example, for a 
neutron star with mass M = IAMq, radius R = 10 km, and the internal temperature T = 10^ K, 
the damping time for the fundamental mode with / = 500 is ~ 160 days. When the temperature 
decreases to T ~ 10^ K, the damping time falls to ~ 15 days. 

In our calculations, we have used the model of the ground-state matter in the neutron star crust 
with a smoothed dependence of the parameters of atomic nuclei on the density of the matter. More 
accurate calculations would require a more detailed model of the ground-state matter, where the 
composition varies with depth in a jump-like manner (showing a series of weak first-order phase 
transitions at some depths). The presence of these jumps could amplify the damping of the 
oscillations. We plan to consider this problem in the future. 
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Appendix. The integration of the viscous dissipation rate of oscillation energy 

over angular variables 

While calculating the angular integral in Eq. (j63|) . it is convenient to introduce the notation 

^ _dU dUk 

c^ik = \- ^— • (68, 

OXk OXi 



Then, the part of Eq. ()63|) . associated with the shear viscosity, can be written as 

dE__l 
'dt ~~4: 



J%r'drJ a,,al,-^\VUf dQ = ^r' (^h - p,^ dr, (69) 



where dfl is a solid angle element, 

Ii = jaik^l- dn and h = j \VUf dQ. (70) 

Here, we have assumed that the unperturbed star is spherically symmetric, so that the shear 
viscosity does not depend on angular variables. The integrals Ji and I2 have been taken ana- 
lytically for the hydrodynamic velocities of the form U = V0, where the velocity potential is 

The integral Ji can be taken if we write the components of the tensor aik in spherical coordinates 
(see, e.g., Ref. ^5). After this, the angular integration is carried out analytically (using the 
properties of the function Yim{0,ip); see, e.g., Ref. |22])- This yields 

J, . 4((F"f + 2i±i(i±il(Ff-6^F'F + ,(, + l)i±i(^ 



where the last equahty is vahd in the approximation of plane-parallel layer. 

Let us now consider the integral 12- For this purpose, we write the divergence of the velocity 

XV^ or or / 



F"+^-l^-^-^F]YU0,v)e-\ (71) 



where An is the angular part of the Laplacian. The integral I2 is easily calculated, 

2 / i/r , 1 \ \ 2 



(72) 



where the last equality is valid in the plane-parallel layer approximation. As expected, Eq. (jHH) 
does not depend on the azimuthal number m (due to the spherical symmetry of the unperturbed 
star). It can be shown that dE/dt > for all allowed values / = 0, 1, 2, . . . 
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